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Optimal control theory is used to develop three automatic 


control systems for the longitudinal motion of a semisubmerged 


ship. A linearized mathematical model of the ship motion is 


used and the control problem is treated as a linear regulator. 


Simulations of the ship's longitudinal motions, utilizing the 


three control systems, are compared for various sea conditions. 


Tt 
Or 
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is concluded that if the wave forces and moments are known 
are estimated as functions of time, a suboptimal controller 
the best controller. However, if a frequency domain approach 
used to estimate the wave forces and moments, the complexity 


a suboptimal controller approaches that of an optimal con- 


troller and the difference between the two is negligible. 


Preliminary work on three techniques for estimating the forces 


and moments is presented. 
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Li" INTRODUCTION 


In 1968 the need arose for a small, relatively inexpensive 
ship which could travel as fast as large Naval vessels and 
provide support for one of the Naval Undersea Research and 
Development Center (NUC) projects -- an unmanned undersea 
vehicle. The ensuing investigationconducted by NUC resulted 
in the semisubmerged ship concept shown in Fig. 1 (Ref. 1]. 

Referring to Fig. 1, the concept is described as follows: 

a. The two streamlined hulls provide the primary 
buoyancy and have less drag and experience smaller wave- 
induced forces than conventional ship hulls. 

E. ‘A propeller 1s placed at the end of each hull 
лесе efficiency is increased by utilizing the boundary layer 
inflow. 

c. Two vertical, surface-piercing bouyant struts are 
attached tozeach hullsspreviding hydrostatiece stability in piten 
and heave due to the waterplane areas and spacing. Hydro- 
dynamic stability in yaw is obtained from the hydrofoil behavior 
of the struts. 

d. The rudders are placed in the propeller slipstreams 
which provides additional rudder force and allows the rudders 
to meee forces on the aft struts. 

e. A horizontal tail surface provides pitch damping as 
well as pitch stability at high speeds. 

AENDE Superstructure, aside from providing the Tormal 


equipment and personnel spaces, provides the juncture fo. the 
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(8а attachment of the two hull-strut systems. This rigid 
system provides stability in roll. 

g- Ine elevator control surfaces incorporated in the 
horizontal tail surface allow a means of controlling the ship's 
response to waves when used in conjunction with 

h. the canard control surfaces which are mounted on 
е forward part of the submerged hulis. 

The purpose of this investigation was to consider various 
ways of controlling the elevator and canard surfaces and make 
recommendations concerning the method of control and areas 


requiring further investigation. 





II. PROBLEM FORMULATION 


А.  NON-DIMENSIONALITY 

The response of a ship to forces -- environmental and/or 
self-induced by control surface deflections -- can be deter- 
mined from a system of equations which describe the ship. If 
the particular ship's displacement is changed, but the hull 
configuration is retained, it is necessary to determine a new 
System of equations in order to determine the ship's response. 
A non-dimensional system of equations however, une allow the 
response to be determined from this one system of equations 
and scaled to any particular displacement. The system of equa- 
tions in this investigation are written in a non-dimensional 


form based on Froude number scaling. 


Froude number = — ІІ, 
vg v 
where 
U = ship's forward velocity 
g = acceleration of gravity 
V = ship's displacement volume 


The dimensional values can be determined by multiplying 
the non-dimensional values by their appropriate factors as 
given in Table I. The Symbol p indicates water density. The 
E с length, le, is the cube root of the displace- 
ment volume, V. The ship's displacement weight is denoted by 


the symbol A. 


О 





TABLE I 
CONVERSION FACTORS 


(Non-dimensional to Dimensional) 


Variable АС а lication 
| kactor 


Length 





Area 

Volume 

Linear Velocity 

Time 

Angular Rate, Frequency 


Linear Acceleration 


Angular Acceleration 


Force 


B. SEMISUBMERGED SHIP FORM USED FOR TREORERTCAL ANALYSTS 

The geometric form used in this discussion and in investiga- 
tions conducted by Naval Undersea Center (NUC), San Dieco, 
California, is shown in Fig. 2. The dimensions are given IN 
terms of the characteristic length, le. Figure 3 will aid in 


interpreting the results for a specific ship size. 


C. COORDINATE SYSTEM 


A body moving in a fluid can move in all six degrees of 


freedom of motion -- translation along each of three ort’ 3gonal 


IT 
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Figure 2.  Semisubmerged Ship Form Used for Theoretical 
Analysis [Ref. 3]. 








Displacement s a long tons 


(Ship length м 5у?) 


Figure 3. Basic Scaling Parameters as a Function of Ship 
Displacement [Ref. 3]. 
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axes and rotation about each of the three axes. The choice 

of a coordinate system should result in one in which the 
equations of motion can be most conveniently analyzed. 

The ship under VUE JOE on has a plane of symmetry (as have 
most ships in existence), which is defined by the ship's 
centerline and a line perpendicular Во ERE ship's deck. Two 
of the axes are chosen to lie in this plane of symmetry with 
the remaining axis being orthogonal to the plane. This choice 
simplifies the expressions for the hydrodynamic forces through 
symmetry and simplifies the equations of motion because the 
axes are parallel to the principal axes of inertia [Ref. 2]. 
If the axes are located such that their origin coincides with 
the ship's center of mass, the axes become the principal axes 
of inertia. 

The ship's position is T in an Earth-fixed right-hand 
cartesian coordinate system. The attitude of the ship can be 
given by an Euler angle transformation from axes parallel to 
the Earth-fixed axis system to the ship-fixed axis system as 
chosen above. Figure 4 shows the relationship between the two 
Systems. 

1: X. - Axis 

The Zain axis is the ship-fixed longitudinal axis in 
the plane of symmetry with the positive direction being toward 


the bow. It is parallel to the ship's keel. 
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У 0 Earth - Fixed 





Ship - Fixed 


Figure 4. Coordinate System 


2. Ya - Axis 
The Yos axis is the ship-fixed transverse axis per- 
pendicular to the plane of symmetry. The positive direction 
is measured to starboard. 
E 2 - Axis 
The Z - axis is the "vertical axis" in the plane of 
symmetry. The positive direction is measured downward toward 
the keel. 
4. Angles of Rotation 
The axis system that has been chosen forms a consistent 
right-handed coordinate system. If $ is the roll angle, 6 15 


the pitch angle and y is the yaw angle, then the positive rota- 


tions are as indicated in Fig. 5. 


D. LONGITUDINAL RESPONSE PROBLEM 
Previous investigations of the semisubmerged ship rest ase 


[Ref. 3] have indicated that the pitch and heave motions in 
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Figure 5. Angles of Rotation 


regular head waves are quite small. However, the ship's 
natural longitudinal response to following waves having a 
height-to-wavelength ratio of 0.05 has been found to be greater 
than desired for wavelengths greater than one ship length. In 
EsEscular, withoaonensdimenszonal speed of 1.65 (high speed), 
the heave response to a wavelength of 10 (height of 0.5) has 
been found to be greater than 1.25. The values for a 2000-ton- 


displacement ship (refer to Fig. 2 for dimensions) would be 


wavelength, Л = 410 ft. 
wave height, h = 20.5 ft. 
heave AZ 51-25 ft. 


This indicates that the ship is thrown out of the water, then 
Slammed back into the water with the upper platform submerging. 
This is an intolerable situation and calls for the design of a 


control system to reduce the ship's natural response. 
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III. SIMULATION MODEL 


A. LINEARIZED EQUATIONS OF MOTION 

The force and moment equations for the ship are written 
in the ship-fixed axis system previously described. If it is 
Е лед that the ship is running at a constant forward veloc- 
ity, U, and that all other motion variables (velocities, angular 
rates and attitude angles) make only small perturbations about 
straight equilibrium running, then the equations of motion can 
be linearized and separated into two groups: (1) the lateral- 
direction equations, and (2) the longitudinal-direction equa- 
tions [Ref. 2]. 

The longitudinal-direction equations of motion [Ref. 3] 
are written below in non-dimensional form. The vertical force 
is given by Eq. (3.1), the pitching moment by Eq. (3.2) and the 
kinematic relationships by Eqs. (3.3). The notation is de- 


fined in Table II. 


(Sala) 
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TABLE II 


DEFINITION OF NOTATION FOR EQUATIONS OF MOTION 


ІТ 527 

Non-dimensional wave forces along the X. axis 
Non-dimensional wave pitching moment about the У axis 
Non-dimensional water plane area 


Non-dimensional X. coordinate of the centroid of water 
plane area 

Non-dimensional metacentric height in pitch 

Non-dimensional radius of gyration about the pitch axis, 
Ys 


Pitch angle 


Non-dimensional position of ship center of mass in 
Earth-fixed coordinates 

Elevator and canard surface deflection angles 

Velocity component along X. and Z, axes respectively 

Angular rate about Y. axis 

Force along 2 axis 

Moment about Y. axis 

Non-dimensional form of uU 

Non-dimensional time, cA 

do/dt 

9 


2 , М 


(1/2 ри у 


ол о> 
О 


7 2/3 
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C 
2, 


од о> 
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2 ^2 ^ ^ 
Ea" = + 
ру4 EUR ICA WC. d 
w q 


+C бе+С ô ] - 2 
mg m. C mp 
= C (3.2) 
2 S 2 En ! 
+А от ү72 + М^ Ww 
" gs ^^ 
а + М, 
9' = а 
и SOT 


The equations of motion can be rewritten in a form where 
they are functions of pitch (9), heave (2), their respecte 
first and second derivatives, control surface deflections and 
wave forces and moments by utilizing Eqs. (3.3). These equa- 
tions are shown below with the coefficients replaced by their 
numerical values as determined by previous investigations 


et. 3]. 


= 1/2 U[-4.303(2'4U0) 


N> 
| 


Са БОЕ 6 

x К (3.4) 
-0.57U 8 1-0.6082 
-0.06470-1.075(2'' 4U0!) 


o ET 
W 


1.4120' = 1/2 O[-4.763(2+U6) 
-15.080'-2,850 8. 
+0.900 6 1-1.0450 | (3.5) 

-0.06472–0. 364 (2'' +00') 


-2 1420" 4M 
wW 
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B. TRANSFER FUNCTIONS 

Reference [3] gives the response of the ship in pitch and 
heave resulting from 0.l-radian control surface deflections. 
These responses can be а from the transfer functions 
relating the responses to the control surface deflections. In 
order to determine if the model used in this investigation was 
equivalent to the one used by NUC, San Diego, the transfer 
functions were derived as shown below and the responses to the 
control surface deflections were obtained. 

The linear Eqs. (3.4) and (3.5) can be Laplace transformed; 
the result of doing this is shown below with the assumption made 
that the initial conditions are zero. The numerical values of 


the coefficients а. апа bi, i=1,...,8, are given in Table III. 


^ ^ ^ 2 
z(s)s а/2(5)5+а,2 (5) +а,8 (5) 5 


+а ,8 (5) га 6 (5) та 6. (3.6) 


За 9 Тава, 


ae b, 8 (s)s+b,6 (s)+b,2(s)s“ 


+b,z2(s)stb.2(s)+b,S, (337) 


+276 +b 


gw 

One method of determining transfer functions is to draw 
A gnal flow graph corresponding to Eqs. (3.6) and (3.7) and 
then apply Mason's Gain Rule to the resulting graph [Ref. 4]. 
Referring to Fig. 6 and Table III for the appropriate cocfficient 


Values, 1t can be seen that the transfer functions wil? ive the 
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6120869 1 b 2 2 U 
-0.2930 b. -0.6701 02-0.2940 
-0.1754 b. 21024 
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+0.4819 b ‚2814 
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Figure 6. Signal Flow Graph. 
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ship's non-dimensional speed, U, as the determining parameter. 


The resulting transfer functions are of the form 


2 
X(s) _ №25 *N,S*NQ i 5 
Y(s) 4 3 2 ` 
D,S +D,S +0 „5 +р 8304 
where the coefficients Ns, i=0,1,2, and Dr 4=0,1,...,4 аге 


given in Table IV. 

Reference [3] investigates the ship's response for two 
distinct speeds, U=1.0 and 1.65. keerring to Fig, 3, itscan 
be seen that for a 2000-ton-displacement ship, these values 
correspond to speeds of 21.8 knots and 36 knots respectively. 
The transfer functions for these two speeds are given in Table 
V, and the responses resulting from 0.l-radian control surface 


deflections were found to agree with those contained in Ref. [3]. 


E COMPUTER MODEL 

The time histories of the ship's responses, pitch, pitch 
velocity, heave and heave velocity, can be found by using a 
computer. A choice, however, must be made whether to use analog 
or digital simulation for the ship model. Both models are given 
below. 

l. Analog Model 

An appropriate model for analog Simulation can be easily 

achieved by normalizing Eq. (3.8) with repect to D4- The result- 
ing model is shown in Fig. 7. A penalty must be paid for the 
ease with which simulation has been achieved and that is that 
the states Xi through X4 have no physical significance. In 
order to determine some physical quantity (e.g., Z, Z', Е, 0"), 
linear combinations of the states Xi through x, are requ. red. 
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Figure 7. Analog Model 


(Assumes no phase reversal through integrators) 


2. Digital Model 

The digital model could easily Be solved by choosing 
the states as being the same as those in the analog model, 
again, however, the chosen states have no physical meaning. 

A more desirable model is achieved by choosing the 
States as being heave, heave velocity, pitch and pitch velccity. 
These states can be obtained by algebraic manipulation of 
Eqs. (3.6) and (3.7). In their present form, each of these 
equations contains the transforms of the second derivative of 
both m and pitch. This explicit interdependence, illustrated 
by the a4-7b4 loop in Fig. 6, must be removed in order to c tain 
a workable model. This is accomplished by substituting i« (336) 
ВЕ (3.7) and Eq. (3.7) into (3.6). 
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The states are chosen as indicated above, resulting 
in a state differential equation which describes the ship. 


This equation has the form 
x(t) = Ax(t)+Bu(t)+Cw(t) (3.9) 


where the state vector is defined as 


7 
ЕЕ 
| 6 
E 
mine control vector is 
| 6, (t) 
u(t) = 
9 (E) 
and the wave-force vector 1s 
Zu 
w(t) = 
M 
ЦИ 
Eho coefficient matrices are 
0 1 0 0 
-0.2951 -0.93620 -0.93620^40.0207 -1.2990 
A = 
Е 0 0 0 1 
0.012 -0.57422 -0.574202-0.2961 -2.0910 
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-0.2964U* 


-0.37060* 


and 


0.4907 
0 


-0.0502 


0 


-0.162402 
0 


0.143202 


305505203 


0.2865 


Equation (3.9) is readily simulated on a digital 


computer to give the ship's time response to control surface 


deflections and known wave forces and moments. 
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MILL CODIUIPRODSSUTSTEMSCEFOR THE SEMISUBMERGED SHIP 


Optimal control theory [Ref. 5] was used to derive three 
controllers for the semisubmerged ship. The derivations of 
the control laws are based upon two assumptions: (1) the 
Ee: in Eq. (3.9) are directly measurable, and (2) no 


measurement noise exists in the measurement of the states. 


A. PERFORMANCE MEASURE 

The performance measure gives a quantitative measure of 
the quality of response of the system under consideration. 
With the exception of long waves (X »14) where the wave crest 
might impact with the upper platform and it would be more 
advantageous to follow the contour of the wave, it is desired 
to maintain the ship as stable as possible in the Earth-fixed 
coordinate system. This requires that the state vector be 
maintained close to the equilibrium point (origin). The con- 
trol surfaces however, should not have excessive deflections 
which would run them into mechanical or electrical Stops, nor 
should they overload the actuation system and result in loss 
of the control system. In view of the above, the requirement 
for the control system becomes: 


The control system must maintain the state 
Wecton, x(t), as close tO the Origin as 
possible without an excessive expenditure 


or control effort. 


The performance measure then chosen is 
E Hi Й 
J = 1/2 J, Ux" (t) gx (t) +" (t)Ru(t) Jat | (4.1) 
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where the final time t; is assumed to be fixed, Q is a real 
symmetric Eve Semsedetinite weighting matrix and R ie a 
real symmetric positive definite weighting matrix. It is 
assumed that the states and controls are not bounded and that 
the final state, x(t.), is free. No penalty has been placed 
on minimizing the deviation of the final state from its 

гей value, 0, It is desired to determine the control 


Et) that minimizes this performance measure. 


B. CONTROL LAW FOR THE SHIP IN A CALM SEA 
The system of equations defining the ship in a calm sea 
(w(t)=0) is 
x(t) = Ax*Bu(t) (4.2) 
The requirement to minimize the performance measure J, given 
by Eq. (4.1) leads to the well-known optimal control law 


[Ref. 5] 


1 


u* (t) -R BK (t) x(t) 


(4.3) 


F (t)x(t) 


И се К(Е 15 the real symmetric positive definite solutionwto 
the Riccati equation 


1 


k(t) = -K(t)A-A'K(t)-Q*K(t) BR B K(t) (4.4) 


which satisfies the boundary condition 


K (te) = 0 


The optimal control law given in Eq. (4.3) is seen to be a time- 
varying feedback of the state vector x(t). The closed-lcop 
System.state equations reduce to 

x(t) = [A+BF(t)]x(t) | (4.5) 
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when the optimal control law is used. Figure 8 shows this 

System in block diagram form. 

fee toe OPTIMAL CONTROL LAW FOR THE SHIP WHEN SUBJECTED TO 
WAVE FORCES 


The Hamiltonian for the forced system is 


(x(t) ,u(t),p(t),t) = 1l/2xl(t)Qx(t)41/2 ul (t) Ru(t) 


+p" (t) Ax (t) +p? (t) Bu (t) 


4p. (t)Cw(t) 
where p(t) is a vector of Lagrange multipliers called costates. 


The theorem that the variation of the augmented functional 


Ј =f 11/2lx (t)0x(t)+u (t)Ru(t) J 
(4.7) 


+p" (t) [Ax (t) +Bu (t) +Cw(t) -x(t) ]}at 


equals zero on an extremal requires that the following neces- 
sary conditions hold for the state equation, costate equation, 


and control [Ref. 51]: 


x*(t) = Ax*(t)+Bu* (t)+Cw(t) (4.8) 
: 9 т fOr all 

p) == ox ^ ox t ADS) t= [torte] (4.9) 
0 = =- За -  Ru*(t)4Blp* (t) (4.10) 


Pere the asterisk denotes an extremal. R is a positive 
definite matrix, therefore Eq. (4.9) can be solved for tke 


optimal control to give 


|: 


u*(t) » -R' Blp*(t) 4.11) 





ISnTeTZuo 5 
э хэварээя Тештз@о чзтм зчета тхозетпбБбән леочтт 





“а етпбтл 
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Mm stitution of Eq. (4.11) into Eq. (4.8) yields 


x*(t) = Ax*(t)-BR™ B p* (t) +Cw(t) (4.12) 


Equations (4.9) and (4.12) now form a system of 2n linear 


nonhomogeneous differential equations 





MEE) См (5) 
----- + |----- 1913) 
р* (+) 0 

хх (Е) Z*E) t. Cw (1) 

----- SO + я ОЕ ст 

p* (t4) p* (t) 0 


(4.14) 


сене ас the state transition matrix of Eq. (4.13). If % 
is partitioned and the integral in Eq. (4.14) is replaced by 
а 2пх1 vector 

te C B 
Па а ат = [----- (4.15) 
t 


Eq. (4.14) can be rewritten as 


129 ИЕ ЕЛЕН ОРЕ (4.16) 


p*(t.) = $51 (grt) g* (0) 4055 (te t) p* (E) FE) (6) (4.17) 


The boundary condition for the costates is 


* 
2 (62) = 0 
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Therefore, Eq. (4.17) can be solved as 


-1 
БЕ СЕ) = m $52 (Е, Е) Фо (Ер, Е) ХР (Е) 
(4.18) 
ML USES 
-22 Е’ -2 
Kalman [Ref. 6] has shown that the required inverse, 
. -1 ў 
922 (te, t), exists for all Е Е [tor tg] SO that Eq. (4.18) 


can be rewritten as 
PAI AE (t) TSt) (4.19) 


Equation (4.19) can now be substituted into Eq. (4.11) resulting 
in the optimal control law (Eq. (4.20)) similar to that in 


Pee. [л]. 


EIS -1 Т 
-R(t)B(t)K(t)x*(t)-R(t)B(t)s(t) 


not) 
(4.20) 


F (t) x* (t) +v (t) 


nt 


Equation (4.20) indicates that the optimal control law contains 
EE Uear function of the system states, and also a command 
сопа], v(t), that not only depends upon the system parameters, 
but also upon future values of the wave forcing function, w(t). 
This indicates that the optimal control law must be anticipatory 
in effect, 

The matrix K(t) and the vector S(t) jean be obtained by 


differentiating Eq. (4.19) 
p(t) = K(t)x*(t)+K(t) x* (t)+s (t) (4.21) 


and substituting Eqs. (4.8), (4.11) and (4.19) into Eq. (4.21). 
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EE = K (t) x* (£) -*K (€) [Ax* (£) +Bu* (t) +Cw (t) ]+5 (t) 


1 


- K(t)x*(t)*K(t) [Ax* (t) -BR ^Bl'p* (t) 4Cw(t) ] 4$ (€) 


!gT (k (€) x* (£) es (£)) 


= K(t)x*(t)-K(t) [Ax* (t) <BR 
+Си (Е) 1+в (+) 
(4.22) 


However, substitution of Eq. (4.19) into the costate Eq. (4.9) 


also results in 


p*(t) - -Qx*(t)-A! (K(t) x* (t) *s(t)) (4.23) 


Equations (4.22) and (4.23) can be equated which, after 


collecting terms, yields 


[K (t) +K (t) A-K (€) BR 1 BTK (€) АТК (Е) +0] х* (Е) 
m (4.24) 
: T = ar > 

II E PERA 5 (Е) (ЕЈ СИЕТ 20 
Ach must be satisfied for all x(t) and w(t). This results 
Meche matrix differential equation for K(t), 

5 T a 

Re HL) AeA KAT) -OTK (EER TB RE) (4.4) 


which is the same as that given earlier for the calm sea 
condition, and an additional vector differential equation, 


S(t) = -[А1-К(&)ВЕ ТВТ] 5 (Е) -К (Е) СИ (Е) (4.25) 


The boundary conditions to be satisfied are 


0 


K(t,) 


S 





Equations (4.4) and (4.25) can now be integrated backward 
in time From te to tg and the resulting solutions can be 
stored and substituted into Eq. (4.20) to give the optimal 
control law 


al 1 


u'(t) = -R'B(t)x*(t)-R B s(t) (4.20) 


= F(t)x*(t)+v(t) 


The controlled system, shown in Fig. 9, is now defined by the 


matrix differential equation 


x(t) = [A*BE(t)]x(t)4Bv(t)4Cw (t) (4.26) 

D. A SUBOPTIMAL CONTROL LAW FOR THE SHIP WHEN SUBJECTED TO 

WAVE FORCES 

The optimal control law given by Eq. (4.20) requires that 
the command signal v(t) be determined either by (1) integrating 
Eq. (4.25) backwards in time assuming the forcing signal w(t) 
iS Клон, or (2) storing previously computed values of v(t) in 
a look-up table and selecting the one which corresponds to the 
appropriate time history of w(t). Either of the above methods 
would require the use of either an existing onboard computer 
Or a computer dedicated specifically to the control system for 
computation and storage. The suboptimal control system .erived 
below was investigated to determine if this added complexity 
could be reduced or altogether avoided. 

Referring to the system equation given by Eq. (4.26), the 
EN ive is to have the state vector, x(t), approach 0. This 
requires that x(t) also approach 0. This can only occur if 


when x(t) O 
Bv(t) = -Cw(t) | (4.27) 


В 
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In general, Eq. (4.27) cannot be solved for v(t) since Bee 


does not exist; however here the first and third rows of 
(Пп В апа С аге zero which allows two new matrices, B' and 


C', to be defined as 


bo) 922 
DAS 
Daz 842 
| (4.28) 
Coi € l 
EE 
а) сб 
It can be shown (see Appendix A) that 
Bi (Е) CE) (4,29) 


EusNSince B' and C* are (2x2) matrices and а exists, Eq. 


(4.29) can be solved to yield 


v(t) o = -B'"e'w(t) 
(4.30) 
= Gw(t) 
The control law can now be written as 
u(t) = E(t)x(t)+Gw(t) (4.31) 
and Eq. (3.9) becomes 
x(t) = [A+BF(t)]x(t)+[BG+C]w(t) (4.32) 


It can be easily shown that as defined, the matrix equality 


BG = -C (За) 
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holds (see Appendix B) so that if w(t) is exactly known, 
Eq. (4.52) reduces to Eq. (4.5), the system equation for calm 


Seas. 
x(t) = [A+BF(t)]x(t) (4.5) 


Therefore, if the suboptimal control law given by Eq. (4.31) 
is used, the computer capability required by the optimal 
control law is avoided. 

Reference [7] shows that the eigenvalues of the matrix 
ENPECC)] have negative real parts so that after initial 
transients, the system states will approach zero. Figure 10 


shows the block diagram for this systen. 


E. CONSTANT STATE FEEDBACK 

The control laws given by Eqs. (4.3), (4.20) and (4.31) 
pintan the time-varying feedback matrix K(t). This matrix 
was determined, as previously stated, by integrating Eq. (4.4) 
backwards in. time. Although K(t) can be computed and stored 
prior to utilization of the above control laws (since the solu- 
tion of Eq. (4.4) does not depend upon the system states), 
the complexity of the control system would be increased cver 
one which incorporated constant rather than time-varying 
feedback of the system states. 

Kalman [Ref. 6] has shown that if (1) the plant is completely 
controllable, (2) there is no cost penalty for final state vector 
deviations from the desired final state, and (3) A, B, R and Q 
Ecuconstant matrices, then K(t) approaches K (a constant 


matrix) as te approaches infinity. To determine whethe: the 
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Plant 1S completely controllable, the rank of the matrix 


ВНЕ AB AB... AB) (4.34) 
must be evaluated. The plant is completely controllable if 
полу if E has a rank of n (in this particular case n is 
four). Equation (4.34) was evaluated and found to have a 
rank of four, so this system is controllable as anticipated. 
The ship meets all of the above requirements with the 
exception that te 


be quite large). However, the solution to Eq. (4.4) using the 


does not approach infinity (although it may 


NPS IBM-360/67 computer (Figs. 11 through 15) shows that the 
feedback gain matrix does reach steady-state constant values 

ша а me short time (= 3.1 sec. for a 2000-ton-displacement 
ШЕР) апа hence K(t) can be approximated by the constant matrix 
[eco Give the following control laws which are used in the 


remaining investigation: 


\ 


пу = кг) (4.35) 
u*(t) - Fx(t) * v(t) (4.36) 
u(t) = Fx(t)+Gw(t) (4.37) 
F can be obtained from 
F = -R nk (4.38) 


where K is the steady-state solution to Eq. (4.4) and hence 
also satisfies the nonlinear algebraic equation 


0 = -KA-À K-QWKBR BK —— (4.39) 
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V. SIMULATION PROCEDURE 


The NPS-360/67 computer was used to solve the ship's time 
response incorporating the control laws determined previously 
for various wave forcing functions. The computer program is 
given following the appendicies, and a flowchart is shown in 
ШЕР ІС: Іп Order to solve Eqs. (4.25) and (4.32), w(t) and 


v(t) must be known. 


A. WAVE MODEL 
The wave forces were estimated using the properties of 


deep water Airey waves. The resulting equations are [Ref. 3] 


2, = Acos(-0'7)+Bsin(-0'7) (5.1) 
М, = Ccos(-o't)+Dsin(-0'7) (5.2) 
where 
с! - non-dimensional wave encounter frequency 
= (2т/Х)1/2 - (2т0/Х) совф (5.3) 
V = wave direction 
Ше (following sea) 
E | (524) 
180° (ahead sea) 


and А is the non-dimensional wavelength of the ocean wave. The 
coefficients A, B, C and D are graphed versus A in Ref. [3]. A 
least squares polynomial curve fitting subroutine, LSQPL2, was 
used to fit a sixth-order polynomial to each of the grap! 3. 


Mach polynomial is a function of М. The resulting polyr nials 
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Figure 16. Computer Flowchart. 
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are used in the computer program to simulate the wave forces 


and moments acting on the ship. 


B. COMMAND SIGNAL 


The optimal control law (Eq. (4.36)) requires the command 
signal 


v(t) = -Е в (6) (5.5) 


Exereo stt) is the solution to Eq. (4.25). Equation (5.5) was 
solved for the three cases shown able VI. All three cases 
are for following seas (y=0%) which, as previously stated, is 
a sea condition ch gives an undesirable ship response. The 
results for case three are shown in Figs. 17 and 18. It is 
seen that as expected for a sinusoidal forcing function, the 


command signal is also sinusoidal in form with the exception 


TABLE VI 


PARAMETERS USED IN SOLVING EQUATION (5.5) 





I= rdentity Matrix 


being where t approaches t, and the boundary condition 
must be met. 
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Since Eq. (4.25 must be solved forward in time and Eq. 
(5.5) can only be solved backward in time, it becomes necessary 
ЕР -пауе a) model төк (Е) to use in solving Eq. (4.25). The 


model that was chosen was of the form 


vi (t) = M,sin(-oT+y,) 
(5.6) 


v(t) = M,sin(-otty,) 

where the magnitudes, Mi and M», and phase angles, Yı and Yor 
were determined from the solutions to Eq. (5.5) and the final 
time deviations from sinusoidal form were ignored. The models 
for the three cases in Table VI are incorporated in the computer 
program as subroutine EST. 

An alternative approach for obtaining the command signal 
is to use a frequency approach instead of integrating Eq. (4.25) 


backward in time. If a new time variable, +, is defined as 


then 
ат = -dt 


and Eq. (4.25) can be rewritten as 
s(t) = [АТ-КВЕ вј5(т) +КСм (т) (5.7) 


which satisfies the initial condition 
s(0) = 0 


ВТК] is stable and since the 


As previously stated, [A-BR- 
transpose of a stable matrix is also stable, Eq. (5.7) is a 


stable equation in t (backwards in time t). 


Sl 





The Laplace transform of s(t), after invoking the 
boundary condition, is given by 


T 


B - 1 
S(s) = [31-АТ+КВВ ВТ)  КСИи (5) (5.8) 


If w(t) is a periodic function with a discrete Fourier spectrum 
Et then the spectrum S{f) is given by 


s -1 
S(£) = L[j2nfI-A -KBR ТВ] KCW(£) (5.9) 


momo it) canbe sOurained by taking the inverse Fourier trans- 
БЕКІП ОС Еа. (5.9). This approach will hold for any general 


Peon) canes the initral™condition of Eq. (5.7) is zero. 
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VI. RESULTS 


- 


The three controllers derived in the previous chapter will 


be called (1) the "constant-gain controller" -- given by Eq. 
(4.35), (2) the "optimal controller" -- given by Eq. (4.36) 
and (3) the "suboptimal controller" -- given by Eq. (4.37). 


Although the constant-gain controller was derived for the 
· сајт sea condition, its use was also investigated for the 
Exundotion where w(t) 0. sure: 19 and 20 show the ship's 
heave and pitch response to a following sea of wavelength 10 


when the 'ship's speed is 1.65 and no control system is used. 


A. CASES INVESTIGATED 
The following assumptions were made in the first five cases 
listed below. 


(1) The control surface deflections were 
instantaneous. 


(2) Ше wavezstorees and moments, w(t), were non-zero 
and exactly known. 


Ве В re Veeror weighring matrix, о, was equal 
БЕШЕНЕ Теле ем matrix, I. 


(4) Contact of the underwater hulls and control 
surfaces with the air-water interface was ignored. 


(5) Initial conditions which were considered severe 


Z(tg) = 0.0 | 
9 (t9) = 0.5 radians 


9 (t9) = 0.0 radians/sec. 


were used (refer to Fig. 2). 
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In the sixth case it was assumed that w(t) was not known 
exactly; the other assumptions indicated above remained in 
force. 

I. Case I 

The reduction of the control weighting matrix R and the 
resultant effect on the ship's response and control surface 
deflections was investigated. 

2. Case II 

The control surface deflections were constrained to 
mo eee andene ereet on the Ship's response was determined. 

3. Case III 

The ship's response using the Constant-Gain Controller 
was determined. 

4. Case IV 

The ship's response using the Optimal Controller was 
determined. 

5. Case V 

The ship's response using the Suboptimal Controller 
was determined. 
6. Case VI 
ШЕСЕ ЕСЕ О рете ЕЩЕ Еа Отит) ОЛ Ее сари 
response was investigated using the Suboptimal Controller. The 


estimate of w(t) was perturbed by +15 percent of its true value. 


B. THE EFFECT OF VARYING CONTROL WEIGHTING MATRIX (CASE I) 
mie valuessofothe weighting matrix R in Eq. (4.1) determine 


how much emphasis is placed upon the expenditure of control 
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effort as compared to maintaining the state vector, weighted 
МС осе со гле origin. It is seen that as the values of 
EXocresuade smaller compared to those of О, more emphasis will 
be placed upon maintaining x(t) close to the origin. 

The ship's response in heave and pitch is shown in Figs. 
„Штала 22 гор а following sea with a wavelength of 10 and a 
ship speed of 1.0. The Optimal Controller was used with two 


Enptereont R matrices 


"Jj 
И 


and 


ye) 
И 


The figures show, as expected, that heave and pitch are reduced 
Ns reduced? The corresponding increase in control effort 
is shown in Figs. 23 and 24 for the Elevator and Canard 


respectively. 


C. CONSTRAINED CONTROL SURFACE DEFLECTIONS (CASE II) 

A large canard deflection (see Fig. 24) was required to 
шоа ties State yector to its Steady-state value when R was 
reduced. The ship's response to the same conditions as in 
Esse aboye with R = 0.52 is shown in Figs. 25 and 26 where 
the control surface deflections were arbitrarily limited to 
un he difference in state trajectories between the un- 
constrained and constrained control surface deflections 3 


almost negligible; additional verification is provided г the 
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performance measure (see Table VII). It was found that an 


additional reduction in the control weighting matrix to 


brought the state vector much closer to the origin but with 
unconstrained control deflections, duced unrealistic deflec- 
tion angles in the initial control period. The control surfaces 
were then constrained to the previous limits and the resultant 
differences in the state trajectories were negligible since 

the period of constraint was small compared to the total control 
interval [tg, tel. As a result, all remaining cases have con- 


ПК ес ен surrtace deflections (t 30°) and 


О. 1 0 


179 
|| 


0 И 


to bring the state vector closer to the origin. 


D. THE CONSTANT-GAIN CONTROLLER (CASE III) 

As previously stated, this controller a derived for 
the case where w(t) # 0. However, since it is based upon feed- 
back of the system states, it will tend to reduce any disturbance 
of the state vector from its equilibrium point. Figures 27 
through 29 show the ship's response and control surface deflec- 
tions for a following sea with a wavelength of 10 when the 
ES speed is 1.65. As expected, the state vector's devia- 
tron from equilibrium has been reduced (refer to Figs. 1” and 
20 for comparison) and the Constant-Gain Controller can е 


used without any knowledge of w(t). 
64 
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Е. THE OPTIMAL CONTROLLER (CASE IV) 

The ship's response using this controller is shown in 
Figs. 30 and 31. The ship's speed was 1.00 and the seas 
were following with a Een of 10. The control surface 
deflections (Fig. 32) are seen to be less than those of the 
Constant-Gain Controller as are the deviations of the state 
vector from equilibrium. This is expected since the Optimal 


Controller minimizes the performance measure J. 


F. THE SUBOPTIMAL CONTROLLER (CASE V) 

Figures 33 through 35 show the ship's response and control 
surface deflections for a following sea of wavelength 10 and a 
ship's speed of 1.65. The response and control surface deflec- 
tions for an ahead sea ( = 180°) of wavelength 10 are shown 
in Figs. 36 through 38. It is seen that, as previously stated, 
the state approaches zero after the initial transients. 

Although the states are maintained at the origin after the 
initial transients, the control effort required is greater than 
that of the Optimal Controller and, as expected, the performance 


measure 1S greater (see Table VII). 


С. THE EFFECT OF INACCURATE ESTIMATION OF WAVE FORCE (CASE VI) 
The response of the ship to imperfect estimates of w(t) was 
investigated for the case where only the estimates of the 
amplitudes were imperfect. It was assumed that the wavelength 
and phase were known exactly. | 
The Suboptimal Controller was used to control the ship 


when the amplitude of w(t) was not exactly known. The с ір was 
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None 


Optimal 
(unconstrained) 


Optimal 
(unconstrained) 


Optimal 
(constrained) 


Constant-Gain 
(constrained) 


Optimal 
(constrained) 


Suboptimal 
(constrained) 


Suboptimal 
(constrained) 


Suboptimal 
(constrained) 


AB: VIT 


PERFORMANCE 


to = 0 sec., 


MEASURES 


te = 80 sec, 


Underestimated |3 





9.328477x10^ 


‚5290623x10° 


2 
- 2 ЗВОВЕЗХЕО 


.5468393х107 


2 


.9093990xI0 
2 
.9204965x10 
2 
.2483721x10 
2 
.9344051x10 


‚3556745x10° 


running in a following sea of wavelength 10 at a speed of 1.65. 


The ship's response in heave and pitch are shown in Figs. 39 


pucEdUEhercecthecamplitude of w(t) was overestimated and undai 


estimated by 15 percent. 


It is seen that the estimate error 


has very little effect on the performance of the Suboptimal 


Controller, 
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VII. PRELIMINARY WORK ON WAVE ESTIMATORS 


The estimators stated below are possible candidates which 
provide an estimate Of w(t) required for the operation of 


both the Optimal and Suboptimal Controllers. 


A. AN INPUT-SIGNAL OBSERVER 
The estimator depicted in Fig. 41 consists of a model of 
the plant being run in parallel with the actual plant. The 


state equation of the model is given by 
e(t) = Ae(t)+Bu(t)+G[x(t)-e(t)] (21) 


and the state equation of the plant is given by Eq. (3.9). Іп 
writing this equation, it has been assumed that all states of 
the system are measured. If the difference between the plant 


мис states (x(t)-e(t)) 1S called c(t); then 


х(+) -е (+) 


S (t) 
- A[x(t)-e(t)]-G[x(t)-e(t) ] tcw(t) (7/923) 
= [A-G]o(t)+Cw(t) 


the solution to Eq. (7.2) 15 


t 


ME o+] e 
0 


[А-6] (t-7) 


c(t) = си (т) ат (7.3) 


so that if the eigenvalues of G are chosen such that @ (A-G)t 


dies out quickly as t increases, then c(t) can be closely 


approximated by 


E 
c(t) = Ј ф(Е–т)Си(т)ат 7.4) 
0 
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u(t) 





Figure 41.  Input-Signal Observer 
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ШІ 9 15 Пе State transition matrix of Eq. (7.2). Referring 
К с апа (5.2), 15 15 seen that Eq. (7.4) can be 
rewritten as 
0 
E 
Ce) = mem a, coswT+a,sinut| dt (755) 
0 


| B4cosuttB,sinut 


The coefficients ау, Aor В. апа В. are constants and can be 


taken outside the integral; thus, Eq. (7.5) can be rewritten in 


thes form 
E t 
о, (Е) = MES, f, (t-t) в1 питат 
(7526) 
E | t 
a (tenes 
(i=l, а) 


The integrals in Eq. (7.6) can be precomputed for a time 
interval (0,t,) and stored in a matrix D(t,) so that Eq. (7.6) 


can be put in the form 


= D(t,) (Fer) 


B2 


If the state differences are now sampled and stored at times 
ti, i-l, ..., k, then an estimate of the coefficients O1, Q5, 
В) апа B, can be obtained by determining the minimum-mean- 


square solution of 
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a EE е = 
сМ Во И у)4®2 (7.8) 
Е By Py 
g 9 | A В. 82| 
given by 
E Se" 
а & |9g(t5) 
^| 2 2 : (7.9) 
um 
в, Jo (ty) | 


where Ø * is the pseudo-inverse of Y. 


В. FAST-FOURIER TRANSFORM ESTIMATION 
An alternative approach is to Laplace transform Eq. (7.2) 


to give | 
0(5) = s [sI-A+G1 е (0) 
E (7.10) 
*[sI-A*G] Cw(s) 


If w(t) is a periodic function, then the spectrum }(f) will be 


-1 
(f£) 7 [j2T7fI-A*G] CW(f) (7 SUL 
шеша вшшасе ог" и(Е) can be determined by 
CW(f) - [j2mfI-A*G]) (£f) (аа 


where )(f) can be determined by taking the Fast-Fourier 
Transform (FFT) of the error between the predicted and actual 


responses. 
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CIN STATE OBSERVER MODEL 
БД Еповс роввірші Бу) 15 to expand the state equation given 
by Eq. (3.9) by defining four additional states. The matrix 
product 
0 
194 cosutta,sinut 
сш) = 0 CS 


а 


1 cosut+Bpsinut | 


can be rewritten as 


0 0 0 0 01 
Cw(t) = coswt sinut 0 0 а. 
0 0 0 0 81 
0 0 coswt sinwt В, | 
(7.14) 
The augmented state vector is now defined as 
z (t) 
z(t) 
Ste) 
x(t) = O(t) (7.15) 
al 
=D 
Bi 
О 


so that the state equation becomes 
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Oo 0 0 0 


| 
| 

e B 
ı coswt -Sinwt 0 0 = 

x(t) = КО 0 0 0 |х(%) %---| ч(%) 
| 0 0 cosut біпше 
СИСЕЕЕ БЕС ГЕЗ = = – – 0 

| D 4 

0 | 0 


(7.16) 


A'(t)x(t)*B'u(t) 
The actual physical states (x4 (t) through x, (t)) are measured. 


The observation matrix for the augmented system is 
H = [I 0] (7219) 


If an observer, Similar to Fig. 41 but with eight states and 
an observation matrix H, is now run in parallel with the 


actual plant so that the observer's state equation is 
e(t) = A'(t)e(t)+B'u(t)+G(t) [He (t) -Hx(t) ] (7.18) 


then the difference between the observer and actual state 


equations can be written as 
o(t) = [A'(t)*G(t) Ho (t) (7.19) 


By choosing G(t) so that [A'(t)+G(t)H] is stable, then the 
SLaues Of therobserver ,~6(t)7 will approach the states x(t) 
NEC LES CE) will approach zero. The states ес (+) 
through eg (t) will then be estimates of Qi, Q5, В, апа В 
respectively. 

Utilizing the above method, an example can be given for 
the first-order case where 


x(t) = -ax(t)+acoswt a»0 
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The states are augmented by defining x, so that 
x (t) ‚ = -аху+ (соѕ t)x, 
хо (Е) = 0 
or 
D NEC EE) 
The observer state equation is 
e(t) = Alt)e(t)+G(t)Hle(t)-x(t)] 


where the observation matrix is 
H = [1 0] 
The matrix 
[-atg, coswt 
ES = 
92: 0 
must now be tested for stability which can be accomplished by 


defining a Lyapunov function [Ref. 8] 


V(t) = 1/2 (o4(t) +05 (t)) (7.20) 


and its time derivative 
v(t) = 0, (t)o, (t) +o, (t)o, (t) (7.7239 
For this example, 


V(t) (-a*g01 (t) * (cosut) o, (t)0, (€) +9,0, (€) 0, (8) 


|| 


(-atg,) 02 (6) + (coswttg,) 0, (t) 5 (t) 


where if the values of gi and 92 are chosen such that 





then Ва. (7.20) 15 positive definite and v(t)«0. The System 
is therefore stable and the estimator states will approach 
the plant states as t approaches infinity. 

If this approach is used, additional work will be required 
ES neure AENA En сап ре selected so that EE e EH is 


stable for the ship. 
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VIII. CONCLUSIONS 


A. THE CONTROLLER 

The easiest controller to implement is the Constant-Gain 
Controller since, as previously stated, it in no way depends 
upon any knowledge (past, present or future) of the wave 
forces and moments. This controller would also be the least 
expensive since the other two controllers are essentially the 
Constant-Gain Controller with additional computed control 
signals. On the other hand, both the Optimal Controller and 
the Suboptimal Controller would be more expensive and complex 
than the Constant-Gain Controller but would maintain the state 
vector closer to the origin. 

The complexity of realizing the command signal depends on 
the nature of the wave forces and how they are estimated. If 
w(t) is sinusoidal in form and is estimated directly in the 
time domain, v(t) can be obtained directly for the Suboptimal 
Controller whereas the Optimal Controller requires either on- 
line or precomputed solutions of Eq. (4.25). In this case, 
the Suboptimal Controller is considered to be less complex and 
„сс ЕЛЕ 1. tf, however, w(t) is some general 
periodic function requiring the FFT estimator or if the fre- 
quency domain approach is used for the single sinusoidal w(t), 
then the difference in complexity between the Optimal and Sub- 
optimal Controllers is negligible and the choice of the "best" 


controller will depend upon how close it is desired to riintain 
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the sedate vector to the origin. In either case, the estimator 
output could be set to zero (if the estimates were bad) and 
both controllers would reduce to the Constant-Gain Controller, 


resulting in still satisfactory performance. 


B. RECOMMENDATIONS FOR FURTHER STUDY 
The assumptions used in this investigation concerning sea 
conditions, control surface deflections and measurements result 
in the recommendations given below. 
l. Sea Conditions 
The constant height-to-wavelength ratio and single- 
waequency wave used in this investigation is not indicative 
of the true sea condition which is a spectrum of different 
height-to-wavelength ratios and frequencies [Ref. 2]. A more 
realistic model for sea conditions should be derived and the 
ship's response to this model investigated. 
2. Control Dynamics 
The E response should be investigated for the 
case where the dynamics of the control surface actuation sys- 
tem are included. This liege icc eR assumed an instanteneous 
control surface deflection capability. 
3. Measurement Noise 
The ideal case of zero measurement noise was assumed 
throughout this investigation. Further investigation is recom- 
mended to include realistic measurement noise. 
4. Wave Estimator 
A wave estimator should be investigated for incl: sion 
in the ship model and the response determined using the stimator 


and measurement noises as noted above. 
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5. Weighting Matrices 
Further investigation is recommended in the determina- 
monto rene values Of the weighting matrices Q and R. The 
possibility of exceeding control surface or structural rates 
and Be accompanying effects on personnel was not included 
in this investigation. If further investigation proves this 
to be a problem area, then new values for the state and control 
weights will have to be determined. 
6. Wave Contouring 
The inclusion of a realistic sea condition model will 
change the point where upper-body contact with the wave occurs, 
previously stated as being wavelengths greater than 14. This 
point of contact should be determined and a controller that 
would allow the ship to follow the wave contour instead of 
driving through the wave should be utilized for this case. 
This would result in a linear tracking problem instead of the 


linear regulator approach used in this investigation. 
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.APPENDIX A 


Here it is shown that the matrix equality 


! = ret 
В Е (4.29) 
holds, where 

5:1 5-2 

Bo 
Par  542 
C21 C22 

on 
C41 C 42 


and from which the command signal, v(t), can be obtained as 
v(t) = -B'e'w(t) 
= = = = (4.30) 


The matrix product 


6,1 Боо | V1 (€) 

Biv(t) = 
Da] 242 Vo (t) 
Рота (Е) tb, v, (t) 


Daz v4 (t) *b, v4 (t) 


and 
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E w (t) 
-C'w (t): 


C (2 w,(t) 


C51W1 (€) *€55w, (t) 


КЕ) +С 


Сати 4272 (5) 


Therefore, since it is given by Eq. (4.27) that 


b,1Vi(t)*b55v5(t) = -(C54w, (t) С ри) (Е)) 


and 


(t) +b - (Саум (t) +C, ow, (t)) 


54171 4272 (5) 


which can be written as 


Бо “1 (Е) + 57 (Е) Соту (6) +С. 2м) (6) 
дедо Са1 1 (Е) +С м 284 
or, 
B'v(t) = - C'w(t) 
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APPENDIX B 


Here it is shown that the matrix equality 


BG - -C i и 


holds, which reduces the system equation to 


E > 


ГА+ВЕ (+) 1х (+), (4.5) 


a Stable system whose state approaches zero as t approaches 


infinity, when the suboptimal controller is used. 


22 


42 


loui —— |) 


BG = B(-B' !c') 
= -вв' le 
Daz  7P22 pi^ NN, 
-b41 Pa сау “2 
(551b4257b55b41) 
pog - 529 
=] [Ae 
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